#include <stdio.h>  
#include <stdlib.h> 
#include <math.h>
#include <string.h>

#define CALLOC calloc

/* Globar variables */
//struct SP_2D_TP * sp_glob;
//struct TP_2D_DATA *  dat_glob;

struct SP_1D /* 1D Spline Structure */ {
	int			order;
	int         n;
	double *	coef;
	double *    knot;
};

struct SP_2D_TP  /* 2D Spline Structure */ {
	int			order;
	int         n1;
	int			n2;
	double **	coef;
	double **   knot;
};

struct SP_3D_TP /* 3D Spline Structure */ {
	int			order;
	int         n1;
	int			n2;
	int			n3;
	double ***	coef;
	double **   knot;
};

struct PWL_2D  /* 2D Piecewise linear spline */ {
	int         n1;
	int			n2;
	double **	coef;
	double *    knot_x;
	double *    knot_y;
};

struct TP_2D_DATA /* Tensor Product Data */{
	int			lenx;
    int			leny;
	double *	x;
	double *	y;
	double **   val;
	double *    weightx;
	double *    weighty;
};

struct L2_1D_DATA{  /* Least Squares Data */
	int			nobs;
	double *	weight;
	double *    xdata;
	double *    ydata;
};

void interv ( double * xt, int lxt, double x, int * left, int * mflag );

double bvalue ( double *t, double *bcoef, int n, int k, double x, int jderiv );

void bsplvb (double * t, int jhigh, int index, double x, 
			 int left, double * biatx );

void bsplvd (double * t, int k, double x, int left, 
			 double * dbiatx, int nderiv );


struct SP_2D_TP * sp_2d_initialize(int order, int n1, int n2,
                                 int  * inform);

struct SP_3D_TP * sp_3d_initialize(int order, int n1, int n2, int n3, int * inform);

struct TP_2D_DATA * data_2d_initialize(int lenx, int leny, int  * inform); 

struct SP_1D * sp_1d_initialize(int order, int n,
                                 int  * inform);

double sp_1d_value(struct SP_1D * sp, double x, int deriv);

double sp_2d_value(struct SP_2D_TP * sp, double x,	double y, int der1, int der2);

double sp_3d_value(struct SP_3D_TP * sp, double x, double y, double z, int der1, int der2, int der3);

double ** get_mat(int n1, int n2);

double *** get_mat_3d(int n1, int n2, int n3);

void sp_1d_set_knots(struct SP_1D * sp, double a, double b);

void sp_2d_set_knots(struct SP_2D_TP * sp, double a, double b, 
				  double c, double d);

void sp_3d_set_knots(struct SP_3D_TP * sp, double a, double b, 
				  double c, double d, double s, double t);

struct SP_2D_TP * sp_2d_initialize(int order, int n1, int n2,
                                 int  * inform);

double ** tp_2d_value(struct SP_2D_TP * sp, double *x, 
				  int lenx, double *y, int leny, int der1, int der2) ;

void free_mat(double ** mat);

void free_mat_3d(double *** mat);

void sp_2d_graph(struct SP_2D_TP * sp, char * file_name,
				 double a, double b, double c, double d,
					  int der1, int der2, int num_graph );

double test_eval(struct SP_2D_TP * sp, int der1 , int der2 );

struct TP_2D_DATA * tp_2d_read(char *file_name, int option, int *inform);

double sp_2d_error(struct SP_2D_TP * sp, struct TP_2D_DATA * dat);

void sp_2d_monotone(struct SP_2D_TP * sp_in, struct TP_2D_DATA * dat_in, 
				 int itmax, int * inform);
void sp_2d_monotone_penalty(struct SP_2D_TP * sp_in, struct TP_2D_DATA * dat_in, 
				 int itmax, int * inform);

void get_max_min(struct TP_2D_DATA *dat, double * max, double * min);

double plus(double x, int k);

void  bchslv ( double **w, int nbands, int nrow, double *b );

void bchfac ( double **w , int nbands, int nrow, double *diag );

void  get_L2_1D_spline( struct L2_1D_DATA * data, struct SP_1D * sp );

int min_int(int i, int j);

struct L2_1D_DATA * data_1d_initialize(int nobs, int  * inform) ;

struct SP_1D * spline_1d_copy(struct SP_1D * sp);

void get_L2_2d_spline( struct TP_2D_DATA * data, struct SP_2D_TP * sp);

int b_center(double *v1, double * v2, double *v3, double * pt, double * u);

struct PWL_2D * pwl_2d_initialize(int n1, int n2, int  * inform);

double pwl_2d_value(struct PWL_2D * sp, double x, double y);


struct SP_2D_TP * sp_2d_controlled_surf(struct TP_2D_DATA * data, int num_pts_per_int,
								int out_x, int out_y, double *** valuep, 
								int * inform);

void free_PWL(struct PWL_2D *pwl);

struct SP_1D * sp_1d_controlled_curve(struct L2_1D_DATA * data, 
						int num_pts_per_int, int * inform);

void free_L2_1D_DATA(struct L2_1D_DATA *p);

void free_SP_1D(struct SP_1D * s);

void free_SP_2D_TP(struct SP_2D_TP * s);

void free_SP_3D_TP(struct SP_3D_TP * s);
